function nmse = nmse(x0, y0, z0, mu, k, N_transient, N_compute)
% 归一化平均同步误差 (Normalized Mean Synchronization Error) 
%
% 输入：
%   x0, y0, z0   - 初始条件
%   mu, k        - MCLM 模型参数
%   N_transient  - 瞬态迭代步数
%   N_compute    - 用于计算 NMSE 的采样迭代步数
%
% 输出：
%   nmse         - 计算得到的 NMSE 值

    %  瞬态迭代
    state = [x0, y0, z0];
    for i = 1:N_transient
        [dx, dy, dz] = mclm(state, mu, k);
        state = [dx, dy, dz];
    end

    % 正式迭代并计算 NMSE
    sum_err = 0;
    for i = 1:N_compute
        [dx, dy, dz] = mclm(state, mu, k);
        state = [dx, dy, dz];
        
        % 计算单步归一化误差
        denom = sqrt(dx^2 + dy^2);
        if denom < 1e-12
            % 若分母过小，可根据需要做特殊处理；这里简单跳过或视为0
            continue;
        end
        
        e = abs(dx - dy) / denom;
        sum_err = sum_err + e;
    end

    nmse = sum_err / N_compute;  % 取平均

end
